Adaptive networks: coevolution of disease and topology 
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Adaptive networks have been recently introduced in the context of disease propagation on complex 
networks. They account for the mutual interaction between the network topology and the states of 
the nodes. Until now, existing models have been analyzed using low complexity analytic formalisms, 
revealing nevertheless some novel dynamical features. However, current methods have failed to 
reproduce with accuracy the simultaneous time evolution of the disease and the underlying network 
topology. In the framework of the adaptive SIS model of Gross et al. [Phys. Rev. Lett. 96, 208701 
(2006)], we introduce an improved compartmental formalism able to handle this coevolutionary 
task successfully. With this approach, we analyze the interplay and outcomes of both dynamical 
elements, process and structure, on adaptive networks featuring different degree distributions at the 
initial stage. 

PACS numbers: 89.75.Fb, 89.75.Hc, 87.10.Ed 
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I. INTRODUCTION 

The vast majority of network-based models of dis- 
ease propagation rely on the paradigm of static net- 
works [3 H] • In this framework, the assumption is made 
that the time scale which characterizes the disease prop- 
agation is much shorter than the time scale with which 
the network structure changes. In contrast to static 
networks, some researchers have investigated the phe- 
nomenon of disease propagation on dynamically evolving 
networks, and have revealed new perspectives on the ef- 
fects of concurrent or casual partnerships [SHS], contact 
mixing [7H5], and demographic changes [TO]. In these 
models, however, the rules which govern the evolution of 
the network are independent of what happens on the net- 
work. Mutual interactions between the network topology 
and the states of the nodes are not taken into account. 

Recently, interest has grown for a new class of networks 
known under the name of adaptive networks |1H I12j . 
They are characterized by the existence of a feedback 
loop between the dynamics on the network and the 
dynamics of the network. Among other applications, 
adaptive networks have been introduced in the study of 
contact processes, such as the study of opinion forma- 
tion [TBH^ and epidemic spreading [25H?7j . In epidemi- 
ological settings, the main idea behind models featuring 
adaptive networks is that individuals may change their 
behavior under the threath of an emerging disease [55]. 
For example, healthy individuals may try to reduce their 
chance of catching the disease by adaptively replacing 
their contacts with infectious individuals by contacts 
with noninfectious ones. This may significantly alter the 
structure of the contact network, thus influencing the way 
the disease will spread. 

Being an emerging field of research, the study of con- 
tact processes on adaptive networks still lacks strong the- 
oretical foundations. Until now, the analytic treatment of 
epidemic models on adaptive networks has been limited 
to low order moment-closure approximations [231 [25H?7] . 



Despite their low complexity, these approaches were able 
to predict novel dynamical features, such as bistability, 
hysteresis and first order transitions. However, their sim- 
ple design does not allow to make accurate predictions 
about the time evolution of the system. An integrated 
analytic formalism able to account for the complete time 
evolution of both dynamical elements, i.e. the spreading 
disease and the evolving network topology, is still lacking. 

In this paper, we present an analytic approach with 
the purpose of filling this important gap. Using for its 
simplicity the epidemic model of Gross et al. [23^ as the 
basic framework, we develop an improved compartmental 
formalism in which nodes are categorized not only by 
their state of infectiousness, but also by the state of their 
neighbors. With this tool, we study the interplay and 
outcomes of disease and topology on adaptive networks 
with various initial configurations. Even if we restrict 
ourselves to one particular model, the approach presented 
here is quite general and could easily be applied to the 
study of other contact processes on adaptive networks. 

The paper is organized as follows. In Secjllj we recall 
the model of Gross et al., and introduce our formalism 
in Sec. |III[ Analytic predictions are compared with the 
results obtained from numerical simulations in Sees. IIVI 
and|V] More precisely, we concentrate on the time evolu- 
tion of the system in Se c. |IV[ while its stationary states 
are investigated in Sec. |V[ Finally, we give further re- 
marks on the endemic stationary state of the system in 
Sec. IVIl and summarize our conclusions in Sec. I VII 



II. SIS DYNAMICS ON ADAPTIVE 
NETWORKS 

We will focus on a simple epidemic model on adaptive 
networks introduced by Gross et al. [551. We consider a 
random dynamical network consisting of a fixed number 
of nodes N and undirected links M = {k)N/2, where 
(k) is defined as the average degree (number of links per 



node) of the network. The nodes of the network rep- 
resent the individuals of a given population, while the 
links stand for potential disease-causing contacts between 
pairs of individuals. Two nodes are said to be neighbors 
if they are joined by a link. Neither can a node be linked 
to itself (no self-loops) nor share more than one link with 
another node (no repeated links). The set of probabil- 
ities {pk{t)} that a node chosen at random at time t is 
of degree fc, called the degree distribution, characterizes 
the topology of the network at this particular time. The 
mean degree of a network corresponds to the first moment 
of its degree distribution, (fc) = ^j, kpk{t) = 2M/N. 

We consider a case of Susceptible-Infectious- 
Susceptible (SIS) dynamics. At any time, each 
node is in a specific state, either susceptible (S) or 
infectious (I). Infectious individuals contaminate their 
susceptible neighbors at rate /3, while they recover 
and become susceptible again at rate a. The coupling 
between disease and topology is implemented by adding 
an adaptive rewiring rule. Susceptible individuals are 
allowed to replace at rate 7 their infectious neighbors 
by individuals chosen at random in the susceptible pop- 
ulation. These rules guarantee that N and M remain 
constant over time. Even if the system contains three 
dynamical parameters, its behavior is characterized by 
two independent dimensionless ratios, e.g. /3/a and 
7/a, since time can always be rescaled according to one 
parameter. 

To perform Monte-Carlo simulations of epidemic prop- 
agation on a network, one requires an explicit knowledge 
of the network structure. Our networks are generated 
according to the following algorithm [2S] • We first gener- 
ate a random degree sequence {ki} of length N subjected 
to the initial degree distribution specified by {pk{0)}. In 
this process, we make sure that ^^ ki is even since each 
link consists of two "stubs" . For each node i, a node with 
ki stubs is produced, then pairs of unconnected stubs are 
randomly chosen and connected together until all uncon- 
nected stubs are exhausted. Afterwards, we test for the 
presence of self-loops and repeated links. All faulty links 
are removed by randomly choosing a pair of connected 
stubs and rewiring them to the former stubs. 

Monte-Carlo simulations of SIS dynamics on adap- 
tive networks are carried out using discrete time steps 
of length At. At each step, the recovery, infection and 
rewiring events are tested with probabilities aAt, PAt, 
and 7At respectively. Self-loops and repeated links are 
explicitly forbidden during the rewiring process. All sim- 
ulations start with the random infection of a fraction e 
of the individuals in the network. We use the parameters 
Ai = 0.1 and N — 25000 in all simulations. The recovery 
rate a — 0.005 is used unless explicitly noted. 

In what follows, we perform simulations on adaptive 
networks featuring different initial degree distributions. 
The first distribution to be used is the delta distribution. 
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each node has the same degree fco. The second type of 
distribution considered is the Poisson distribution. 



Pk 



z^e 



k\ 



(2) 



which corresponds, in the limit N ^ 1, to networks in 
which the presence of a link between two nodes is gov- 
erned by the same probability, independent of the links 
already present in the network. We will refer to them as 
Poisson (P) networks, and their mean degree is given by 
{k)^ = z. Finally, we will also use a truncated power law 
distribution, 
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< fc < fee 
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where t > and C — 2_/feLi ^^^ so that the distribu- 
tion is properly normalized. This produces power law 
distributed (PL) networks, where highly connected hubs 
and individuals with few connections coexist. The mean 
degree of networks generated by Eq. (l3| is given by 

{k)^^ = ^~^ X^feLi ^^~^- To obtain a network with 
(/c)^^ = 2, we use r = 2.16104 and fc^ = 20. 



III. IMPROVED COMPARTMENTAL 
FORMALISM 

In order to describe the complete time evolution of 
the model defined in the last section, we introduce an 
improved compartmental formalism in the spirit of the 
formalism presented in the Appendix of Noel et al. [30] . 



A. Dynamical equations 

Let Ski{t) and Iki{t) be the fractions of nodes of total 
degree k and infectious degree I < k that are respectively 
susceptible and infectious at time t [21] ■ Here by total 
degree we mean the total number of links that belong 
to a node, and by infectious degree the number of those 
links shared with infectious individuals. We define the 
zeroth order moments of the Ski and Iki distributions by 



5 = ^ Ski and I = y^Iki 
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the first order moments by 

55 = 5I(^-0^H > Sl = Y.lSk. 
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Is = ^{k~ l)Iki and h = ^ Hki 
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which produces a degree- regular (DR) network where 



and the second order moments by 
Ssi = ^(fc - l)lSki , Sii = Y, 1(1 - ^)Ski , etc. (6) 
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FIG. 1. Schematic illustration of the events described in the 
text which result in changing a node from one compartment to 
another. Susceptible nodes are represented by open symbols 
(O) and infectious nodes by filled symbols (•). 



Physically, the zeroth order moments correspond to the 
density of S and I nodes, the first order moments to the 
density per node of the various types of arcs, and the 
second order moments to the density per node of the 
various types of triplets in the network [32'. 

As mentionned in the last section, our model is con- 
strained by two conservation relations, namely the con- 
servation of nodes. 



S+I=l , 
and the conservation of links, 

Ss + Si + Is + Ii = (k) , 



(7) 



(8) 



which must hold at any time t. Moreover, since the net- 
work under consideration is undirected, the density of 
SI links must always equal the density of IS links. This 
yields the additional constraint 



Si 



(9) 



We now derive an ordinary differential equation (ODE) 
for each compartment of the system. Let us illustrate the 
reasoning for Ski- First, nodes can change compartment 
according to a change in their own state. Nodes are added 
to Ski at rate ahi as nodes from hi recover and become 
l(a)|, and are removed from Sm 

-^ -~ T(b)]. 



susceptible again [Fig 

at rate /3lSki as they get infected [Fig 

Second, we have to account for a change in the state 
of a node's neighbors. Nodes from Ski are transfered to 
Sk{i-i) when one of their infectious neighbors becomes 
susceptible again [Fig. 1(c) , which occurs at rate alSki- 
On the opposite, nodes from Ski are moved to Sk(i+i) 
when one of their susceptible neighbors catches the dis- 
ease [Fig. l(d)| . Since we do not know the exact number 
of infectious neighbors connected to each of the (fc — I) 
susceptible neighbors of nodes in the Ski compartment, 
our best option is to take an average over the entire 



network. In doing so, we make the assumption of zero 
degree correlation: we assume that nodes, beyond the 
knowledge of the state of infection of each neighbor, are 
connected at random in the network. Hence, a node in 
compartment Sk'i' with V infectious neighbors will be 
reached by following a SS link with a probability equal to 
{k' — l')Sk'i' I Ss- Thus, Ssi/Ss gives the average number 
of infectious neighbors that a susceptible node reached 
by following a SS link has. The rate associated to the 
Ski -^ Sk{i+i) transition is therefore P{Ssi/Ss){k-l)Ski. 
Finally, we have to account for the effects of rewiring. 
Nodes from Ski become labeled as Sku-i) at rate jlSki 
as they break connections with their infectious neighbors 
[Fig. 1(e) . Moreover, a node from Ski is moved to the 



compartment S^ 



ik+l)l 



if it is chosen as the "new neighbor" 



in a rewiring event [Fig. l(f)|. Since the strength of 



rewiring events is jSj and a node from Ski is randomly 
chosen with a probability Ski/S, this occurs at a rate 
^{Si/S)Ski- By summing all contributions, we obtain 
the following ODE governing the time evolution of the 
Ski compartment: 
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dSk. 
dt 



= ahi - PlSki +a {1 + l)Sk{i+i) - IS, 
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{k~~l + l)Sk{i-i) - (fc - l)Ski 
Si 
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A similar reasoning for the hi compartment yields the 
following ODE: 



dhi 
dt 



ahi + PlSki + a 
Sii 



+ /? 1 + 



(/ + l)4((+i) - Ihi 
{k^l + l)h(i-i) - (fc - l)hi 



■7 



{k-l + l)I(k+i)i -{k- l)Ii 
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(11) 



It is straightforward to show that the infinite system of 
ODEs consisting of Eqs. (llOl) and (11) satisfies the con- 
straints given by Eqs. (|7|-19|. 



B. Initial conditions 

In order for the dynamics to be completely specified, 
we need to write an initial condition for each compart- 
ment. In the case where a fraction e of the nodes is 
initially infected at random, they are given by 

Sk0) = {l-e)pk{Q)i^^e\l-ef-' (12) 

and 

hm^epk{0){^^e\l-ef-' . (13) 

Again, we easily verify that this set satisfies Eqs. (|7])-(|9|. 
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FIG. 2. (Color online) Disease prevalence 7 against time t on networks featuring the same mean degree (fc) = 2 but different 
initial degree distributions {p^^, p^ and pf^). Parameters are (a) P — 0.06, 7 = 0.02 and (b) /3 — 0.04, 7 — 0.04. We use 
a = 0.005 and e = 0.1 in all simulations. Points and error bars (if larger than marker size) correspond to the mean and standard 
deviation computed over 10000 Monte-Carlo simulations; solid lines are the predictions of our analytic approach. 



The complete time evolution of the system, is obtained 



by integrating numerically the set of ODEs given by ( 10 ) 
and (111 truncated at fcmax, together with the initial con- 
ditions (12) and (13 1. Constraints (|7])-(I9]) can be used to 
check the precision of the numerical integration. The 
complexity of the system of equations is 0{k^g^^). 



IV. TIME EVOLUTION: INTERPLAY 
BETWEEN DISEASE AND TOPOLOGY 

As stated previously, we initialize the dynamics of the 
model by infecting a fraction e of the nodes in the network 
at random. Afterwards, the states of the nodes and the 
topology of the network coevolve according to the rules 
prescribed in Sec. In] In this section, we analyze the time 
evolution of the system, both from the perspective of the 
spreading disease and of the evolving network topology. 

A first quantity of interest is the evolution of the dis- 
ease prevalence, defined as the fraction of infectious indi- 
viduals at time t. According to our previous definitions, 
it is simply given by /. In Fig. [2] we illustrate the evo- 
lution of / for networks with different initial topologies, 
namely DR, P and PL networks featuring a mean de- 
{k) 



gree 



2. In Fig. 2(a 



all systems reach an endemic 
steady state where the disease prevalence seems to sta- 
bilize at the same value. However, the time evolution of 
/ follows very different patterns depending on the initial 
configuration of the system. In Fig. 2(b) we show that 
for another set of parameters, the P and the PL network 
converge towards an endemic state, whereas the rewiring 
is sufficiently strong to hinder the initial propagation in 
the DR network and consequently, the system converges 
to a disease-free state. This indicates that the initial net- 
work topology may also influence the global outcome of 
a particular epidemic scenario. 

As the disease propagates, connections between indi- 
viduals are being adaptively rewired, which affects the 



degree distribution of the network. The normalized de- 
gree distributions of susceptible and infectious individ- 
uals are given in our formalism by s^ = X); Ski/S and 

We consider the simplest example of a population ini- 
tially connected via a DR net work w ith (fc) — 2. For 
the same parameters as in Fig. 2(a) we show in Fig. ^ 
the time evolution of the probabilities Sk and ik for low- 
degree (fc = 0, 1, 2) and high-degree {k = 3, 4, 5) nodes. 
At i = 0, all nodes are of degree 2. When the dis- 
ease starts to propagate, both degree distributions are 
rapidly modified. The fraction of degree 1 infectious 
nodes quickly increases, because the susceptible neighbor 
of a degree 2 node who has just been infected will try to 
rewire its connection. This result in an increase of degree 
3 and higher susceptible nodes. During this first phase of 
infection, we observe in Fig. |3(b)] that the fraction of high 
degree susceptible nodes smoothly increases, each degree 
probability lagging behind the preceding one. Shortly 
after t = 200, the fraction of susceptibles that are of de- 
gree 3 and higher suddenly drops, while the fraction of 
them that are either of degree or 1 increases. Since sus- 
ceptible nodes cannot lose connections, this means that 
the infection reaches an important number of high de- 
gree susceptible nodes that have formed during the ini- 
tial stage of infection. Afterwards, the topology settles 
slowly towards its stationary state. 

In order to investigate further the interplay between 
disease and topology, it is useful to look at the evolution 
of some important observables of the system. First, we 
consider the fraction of SI links in the network. This 
quantity, given by Sj (or Is), is directly proportional to 
the number of new infections at a given moment in time. 
While Si is a good measure of the instantaneous danger- 
ousness of the situation, it does not yield any information 
about the potential of the disease to spread further. To 
quantify the latter effect, we use the effective branching 
factor Kjg = Ssi / Si, which is the average number of sus- 
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FIG. 3. (Color online) Degree probability s^ and ik for susceptible (open symbols, dashed lines) and infectious (filled symbols, 
solid lines) individuals against time t on an adaptive networks with a DR initial degree distribution (fco = 2). The parameters 
of the system are a = 0.005, j3 = 0.06, 7 = 0.02, and e — 0.1, as in Fig. |2(a)[ Points and error bars (if larger than marker size) 
correspond to the mean and standard deviation computed over 10000 Monte-Carlo simulations; curves are the predictions of 
our analytic approach. 



ceptible neighbors of a susceptible individual reached by 
following an IS link. Since n^g is related to the degree 
of correlation between susceptible individuals, we also 
compute the average fraction Css = Ss/{Ss + Si) of 
connections that susceptible individuals share with other 
susceptible individuals. 

In Fig. Wl the time evolution of S"/, K^g and Css is il- 
lustrated lor the example of the initial DR network con- 
sidered so far. We now put these results in parallel with 
those obtained on Figs. |2(a~ andls] In the system under 
consideration, /3 and 7 arc significantly larger than a, and 
therefore, the infection and rewiring processes dominate 
at early times. At f = 0, everyone is of degree 2, which 
means that infecting a susceptible node does not increase 
the number of SI links. Hence Sj quickly decreases due to 
adaptive link rewiring and the propagation speed gradu- 
ally decreases. As we saw in Fig. l3J initial rewiring events 
result in an increase of high-degree susceptible nodes. 
Morever, we see in Fig. |4] that those susceptible nodes 
form a strongly linked community: near t = 100, Css 
stays as high as 0.9 even while nearly half of the nodes in 
the network are infected. The initial disease propagation 
phase is thus characterized by a topological segregation 
between susceptible and infectious individuals. 

However, this situation is unstable. There still are 
some SI links in the network, and a strongly connected 
community of susceptible nodes is very vulnerable to a 
rapid infection. The disease eventually invades the sus- 
ceptible community: near t = 100, the effective branch- 
ing factor Kjg rises again, yielding an increase in the 
number of SI links which reaches a peak near t = 200. As 
the disease propagates in the community of susceptible 
nodes, Css exhibits a sharp decrease. To this invasion of 
the tightly linked community of susceptible nodes corre- 



spond the second burst of infection observed in Fig. 2(a) 



and the sudden decrease of high-degree susceptibles ob- 
served in Fig. [3J both near t = 200. After this sec- 
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FIG. 4. (Color online) Time evolution of the fraction of SI 
links Si, the effective branching factor kJs, and the aver- 
age number of connections that susceptible nodes share with 
other susceptibles Css in a system with an adaptive networks 
featuring a DR initial degree distribution (fco = 2). The pa- 
rameters of the system are a = 0.005, /3 = 0.06, 7 = 0.02, and 
e = 0.1. Points and error bars (if larger than marker size) cor- 
respond to the mean and standard deviation computed over 
10000 Monte-Carlo simulations; solid lines are the predictions 
of our analytic approach. 



ond phase of disease propagation, the system converges 
smoothly towards a stationary state. 

Figures [2| [3) and |4] confirm that our formalism is 
well capable of tracking the time evolution of disease 
and topology on adaptive networks. Numerical results 



obtained from Monte-Carlo simulations are in excellent 
agreement with analytic predictions. We may mention 
that a good agreement between theory and simulations 
has also been obtained for the time evolution of degree 
distributions and topological observables in systems fea- 
turing P and PL initial networks with the parameters 
used in Fig. |2(a)| 



V. STATIONARY STATES 

After studying the time evolution of the system, we 
now investigate its stationary states. At first glance, a 
given epidemic scenario may admit three different out- 
comes. First, the disease may not be virulent enough to 
propagate throughout the network, hence the system will 
converge towards a disease-free state, i.e. a frozen con- 
figuration where all the nodes are susceptible. Second, 
the disease may reach and maintain a fixed macroscopic 
prevalence in the population, where the number of new 
infections equals the number of recoveries at any time. 
When this endemic state is reached, the system is in ac- 
tive equilibrium, since infection, recovery and rewiring 
events continuously occur. Third, we cannot reject the 
possibility that the disease prevalence may never settle 
to a constant value, and behave in a periodic, quasiperi- 
odic or even chaotic fashion. In our study, we only re- 
ported the first two scenarios, what is consistent with 
previous studies of this or similar systems [231 I25ti27j . 
The presence of a stable limit cycle in a narrow region of 
parameter space was theoretically predicted in [23], and 
short-lived oscillations were reported in Monte-Carlo sim- 
ulations of large systems ,33] . In what follows, we do not 
report any oscillations. We believe that oscillations could 
in principle also be found with our analytic framework. 
However, due to the high complexity of our approach, 
finding the oscillatory regime would require a systematic 
investigation of the parameter space beyond the scope of 
this paper. 



A. Bifurcation structure and topology at 
equilibrium 

In order to study the properties of the system at equi- 
librium, we first consider the stationary disease preva- 
lence /* . In our analytic formalism, /* is obtained by in- 
tegrating Eqs. ( 10 ) and ( [TTj ) until convergence is reached 
towards a stable manifold. Predictions of our analytic 
formalism are compared in Fig. [5] with the outcome of 
Monte-Carlo simulations for static and adaptive networks 
with an initial Poisson degree distribution. For compari- 
son, we also illustrate the analytic predictions of the low 
order approach of Gross et al. (see Appendix) . 

where each node has a mean degree of 



the correct equilibrium behavior of the system with and 



without rewiring. On Figs. 5(b) and 5(c) we decrease the 



Figure 5(a) 



mean degree of the system to (fc) = 7 and 2 respectively. 
We see that our formalism continues to remain valid as 
{k) diminishes, while the formalism of Gross et al. looses 
its accuracy. This result can be mainly explained by the 
fact that unlike in our approach, the equations of Gross et 
al. do not distinguish between the individual behaviors 
of nodes with different degrees. These behaviors become 
increasingly heter ogeneous as (fc) is decreased [34]. Fi- 
nally, in Fig. 5(c)[ we see that analytic predictions for /* 
are more accurate on adaptive than on static networks. 
This is due to the fact that link rewiring induces a certain 
amount of shuffling in the network connections. In this 
case the history of the transmission events that did or did 
not happen becomes less important and the description 
of the system at a coarse-grained level is more accurate. 

For the systems with link rewiring shown in Fig. [5l we 
can clearly see the existence of a bistable regime char- 
acterized by two first order transitions. To these dis- 
continuous transitions correspond two thresholds: the 
persistence threshold /3pcrj from which an already well- 
established epidemic can persist in the population, and 
the invasion threshold /^inv, where the disease- free state 
become unstable for all finite values of e. These features 
have already been recognized as generic features of epi- 
demic models on adaptive networks [33 US] |3S] . 

After illustrating the effect that link rewiring has on 
the stationary disease prevalence in the system, we can 
also study how it affects the topology of the underlying 
network. The case of higher interest is the topology of the 
endemic state, where the system is in active equilibrium. 

In Fig. [6| we illustrate the normalized stationary de- 
gree distributions observed at various rewiring rates in 
the endemic state of a system with an initial P net- 
work and (fc) = 7. In our formalism, they are given by 
si = ^i Sli/S* for susceptible nodes and i^ = ^; IlJl* 



for infectious nodes 
Fig.g 



For a static network, shown in 
both stationary degree distributions follow a 



Poisson distribution. The peak of i^ is found at higher 
degree than the peak of sj because high degree nodes 
are more likely to get infected. For adaptive networks, 
shown on Figs. |6(b")] and 6(c) both stationary degree dis- 
tributions get significantl y bro ader, particularly for sus- 



ceptible nodes. In Fig. 6(d) the variance a^ of both 



{k) = 20, corresponds to the case treated previously by 
Gross et al. in |23| . In this highly connected limit, we 
see that both analytic formalisms are able to reproduce 



distributions is plotted versus the rewiring rate 7. a 
is a smoothly increasing function of 7 for both distri- 
butions, its increasing rate being greater for susceptible 
nodes. For comparison, we also indicated the variance of 
the stationary distributions on a static network (7 = 0). 
Our results show no apparent continuous transition be- 
tween the degree distributions at equilibrium on static 
and adaptive networks. Starting from the equilibrium 
topology in the adaptive regime, it is therefore impossi- 
ble to recover the initial topology by slowly decreasing 
the rewiring rate of the system. We will return to the 
implication of this observation in more details in Sec. |VI| 
The results presented in Fig. [6] confirm that link 
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FIG. 5. (Color Online) Bifurcation diagrams of the stationary disease prevalence /* versus infection rate /3 on static and 
adaptive networks with an initial Poisson degree distribution and (a) (fc) = 20 and a = 0.002; (b) (fc) = 7 and a = 0.005; 
(c) (fc) = 2 and a = 0.005. Solid lines are the predictions of our analytic approach while dashed curves are the predictions of 
the formalism of Gross et al. (stable and unstable manifolds). Points and error bars (if larger than marker size) represent the 
outcome of Monte-Carlo simulations. To collect statistics in Monte-Carlo simulations, at least 25 simulations were run for each 
value e — 0.001, 0.01, 0.05, 0.99 and 0.999. For each run, the initial transient was discarded, and the prevalence at equilibrium 
was averaged over at least 10000 time steps. Averages were calculated with the sets of results obtained from this procedure. 



rewiring has two main effects of opposite epidemiologi- 
cal consequences [23]. On the one hand, it locally pro- 
motes the isolation of infectious individuals, but on the 
other hand, it triggers the formation of highly connected 
individuals, which acts as superspreaders of the disease. 
This dual effect may be responsible for the apparition of 
a bistable regime in parameter space, which is not ob- 
served in static networks. 

The effects of link rewiring on the topology of adap- 
tive networks have been previously observed in stochastic 
simulations [23, 26^, but until now, no analytic approach 
was able to model them correctly. A previous attempt 
by Shaw and Schwartz |26l has been unsucessful, basi- 
cally because their formalism was not able to account for 
the correlations between the state of neighboring nodes. 
The results shown in Fig. [6] highlight that our improved 
compartmental formalism is able to capture with great 
accuracy the degree distributions of the system at equi- 
librium. By characterizing nodes by their total and in- 
fectious degree, we are able to overcome the correlation 
problems faced in [55] . 



B. Comparison of phase diagrams for different 
initial networks 



In the last subsection, we have studied the behavior of 
the model on adaptive networks with an initial Poisson 
degree distribution. For these particular initial networks, 
we have confirmed that there exists a bistable region at 
finite rewiring rate in parameter space. In this region, 
the initial disease prevalence plays an important role in 
determining if an epidemic will either die out or persist 



in the population. Moreover, we saw in Sec. IV that the 
initial network topology also influences the evolution and 
outcome of a system. We now study the location of the 
persistence and invasion thresholds in systems featuring 
different initial topologies. 

Phase diagrams in the plane (7, /3) for three systems 
featuring different degree distributions at the initial stage 
with (fc) = 2 are illustrated in Fig. [t] The location of 
both thresholds /3per and jSi^v were obtained with our 
improved compartmental formalism using a bisecting al- 
gorithm, with e = 0.0001 for Pi^v and e = 0.99 for /3pcr. 
Figure [7] shows that all three networks display a bista- 
bility region between regions where only one stationary 
state, either endemic or disease-free, is stable. At fixed 
recovery rate a, the extent of this bistability region de- 
pends on the rewiring rate 7 and the initial topology of 
the network. The invasion threshold /3inv grows much 
faster as 7 is increased in systems with an initial DR 
network than in systems with initial P and PL networks. 
These results suggest that for the same link density (fc) 
and at small e, link rewiring as a disease control strategy 
is more efficient in homogeneous networks, i.e. networks 
with small fluctuations in their degree distribution. As 
we mentioned previously, adaptive rewiring tend to su- 
press disease propagation on a local scale, but has the 
potential to create high degree susceptible nodes on a 
global scale, which favors the spreading of the disease. 
On static networks, the initial spreading phase is known 
to be slower on homogeneous networks than on strongly 
heterogeneous networks |35j . When an adaptive rewiring 
rule is added, it is then easier for homogeneous networks 
to hinder the initial propagation of the disease on a local 
scale before it reaches a macroscopic prevalence and a 
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FIG. 6. (Color online) Normalized stationary degree distri- 
butions Sfc and il in the endemic state for susceptible (open 
circles) and infectious (filled diamonds) nodes for a = 0.005, 
13 = 0.008, and (a) 7 = 0, (b) 7 = 0.025, and (c) 7 = 0.050. 
(d) Variance a^ of the stationary degree distributions of sus- 
ceptible and infectious nodes versus the rewiring rate 7. The 
initial network used has a Poisson degree distribution with 
(k) = 7. Points and error bars (if larger than marker size) 
correspond to the mean and standard deviation computed 
over (a-c) 1000 and (d) 250 Monte-Carlo simulations; solid 
lines are the predictions of our analytic approach. 



critical concentration of high degree nodes is attained. 
Consequently, /3inv is higher in homogeneous adaptive 
networks. 

Except at very small rewiring rates, Fig. [7] shows that 
the persistence threshold is the same for the three sys- 
tems. Since this threshold marks the point from which 
a stable endemic state appears in the system, this sup- 
ports the existence of a universal endemic state common 
to those systems, regardless of the initial topology. We 
will return to this point in Sec. |VI| We believe that the 
persistence thresholds differ when 7 is small because the 
network does not evolve rapidly enough in all systems to 
converge towards a topology on which an endemic state 
would be stable, even if such a topology exists. 

The results presented in this section illustrate the im- 
portance of initial conditions in determining the global 
outcome, either endemic or disease-free, of a given epi- 
demic scenario. The initial network topology determines 
the size of the bistability region, and inside the latter, 
the initial disease prevalence determines which station- 
ary state will be reached. 
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FIG. 7. (Color online) Location of the persistence /3por and in- 
vasion /3inv thresholds versus the rewiring rate 7 as computed 
from our analytic approach at a = 0.005 for systems featur- 
ing the same mean degree (fc) = 2 but different initial degree 
distributions {pk^ , Pk a^nd Pt^)- For each initial degree dis- 
tribution, the upper-most region corresponds to the endemic 
state, the lower-most region to the disease-free state, and the 
inner region to the bistable regime. Except for very small 
values of 7, the persistence threshold of all three networks is 
the same (lowest curve). 



VI. FURTHER REMARKS ON THE ENDEMIC 
STATE 

In the last sections, we have gathered much evidence 
supporting the claim that the endemic state found in sys- 
tems featuring an adaptive network is only determined by 
the dynamical parameters of the system, {a,/?, 7}, and 
the link density (fc). For a given set of these parameters, 
the endemic state appears to be universal, i.e. it does not 
depend on the particular initial conditions of the system. 

Let us briefly recall the results that corroborate this 
idea. For systems with the same value of {a,/3,7} and 
(fc) , we showed in Sec. IV that even if their evolution to- 



wards the endemic state, if ever reached, is different, the 
disease prevalence / converges in all cases to the same 
value. In addition. Fig. [8] displays the fact that their 
degree distributions also converge to the same distribu- 
tion. In Sec.lV Al we found no continuous transition in P 



networks between the stationary degree distributions of 
susceptible and infectious individuals as link rewiring is 
turned on. Moreover, we computed in Sec. [VB] the value 
of the persistence treshold /3pcr, and found that it has the 
same value in systems featuring adaptive networks with 
different initial topologies but the same mean degree (fc). 
Is it reasonable to ask if this claim makes sense. In 
the coevolutionary model studied here, susceptible indi- 
viduals are allowed to avoid contact with infectious indi- 
viduals by changing acquaintances. As the disease prop- 
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FIG. 8. (Color online) Normalized stationary degree distri- 
butions si and ife in the endemic state for susceptible and in- 
fectious nodes in three systems with different initial network 
topologies characterized by (k) = 2. Simulation parameters 
are a = 0.005, /3 = 0.006, and 7 — 0.02. Points correspond to 
the mean computed over 1000 Monte-Carlo simulations; solid 
lines are the predictions of our analytic approach. 



agates, links are being rewired and the network slowly 
loses memory of its initial structure. As already men- 
tioned, the endemic state of the system is active. For a 
system able to reach it, the coevolution process between 
state and topology lasts for a infinitely long time. Hence, 
at some point, the information about its initial structure 
is completely lost. 

This phenomenon can be interpreted in the frame- 
work of statistical mechanics. The endemic state can 
be thought as the state of maximum entropy of the sys- 
tem. It only depends on the density parameter of the 
system, (fc), and the interaction parameters between par- 
tides (nodes) of the system, {a,/3,7}. As the system 
evolves dynamically towards the state of maximum en- 
tropy, information is lost. Therefore, the evolution pro- 
cess towards the endemic state is irreversible. 

However, even if there is much evidence in favor of the 
existence of a universal endemic state for given {a, /3, 7} 
and (fc), this statement still remains at the conjecture 
level. Since our improved compartmental formalism does 
not seem to admit an analytic solution for the equilibria 
of the system, it is impossible to mathematically demon- 
strate that the solution for the endemic state does not 
depend on the particular initial conditions of the system. 
Another approach may be needed to solve this problem. 



VII. CONCLUSION 

In the spirit of the formalism presented in Appendix 
of |30j . we have introduced an improved compartmental 
formalism in the framework of a simple SIS model on 
networks featuring an adaptive rewiring rule '23] . In our 



approach, individuals are put in compartments accord- 
ing to their state of infectiousness, their total degree k 
and their infectious degree I. With these considerations, 
a set of ODEs describing the dynamics of the system is 
obtained and can be integrated numerically to yield its 
evolution and stationary states. Theoretical predictions 
were found to be in excellent agreement with numeri- 
cal results for adaptive networks with various degrees of 
heterogeneity at the initial stage. Being the first capa- 
ble of reproducing the complete time evolution of both 
dynamical elements, process and structure, the approach 
presented in this paper marks an important step forward 
in understanding the complex behavior of adaptive net- 
works. 

As a pedagogical example, we have analyzed in details 
the coevolution of disease and topology in a system fea- 
turing a degree-regular network at the initial stage. By 
tracking the evolution of meaningful observables, we were 
able to point out the dual effect of link rewiring in the 
population. Beside bringing better insights about the 
interplay between disease and topology on adaptive net- 
works, this simple example showed that our formalism is 
very well-suited for the study of these complex systems. 

Moreover, the results obtained show that the initial 
conditions, i.e. disease prevalence and network topology, 
play an important role in determining the evolution and 
outcome of a particular epidemic scenario on an adap- 
tive network. It does not only affect the speed at which 
stationarity is reached, but can also determine which sta- 
tionary state is reached - either endemic or disease-free. 
In contrast, the properties of the endemic state do not 
seem to be affected by the initial topology of the network. 
We have presented strong evidence that it only depends 
on the dynamical parameters of the system, {q;,/3,7}, 
and the mean degree of the network, (fc). However, since 
our model cannot be solved analytically, this conjecture 
remains to be proved. 

The use of the model of Gross et al. as the framework 
of this paper served the purpose of proof of concept. De- 
spite its appealing simplicity, this particular model is ac- 
tually lacking in realism. Possible directions for further 
research could consist of including more realistic features 
in epidemic models on adaptive networks. A first step 
in this direction could be the introduction of cliques [35] 
in the network to account for the community structure 
observed in many real- world networks [37] . Another in- 
teresting effect to implement would be some mechanism 
of preferential rewiring. For example, this could be mod- 
eled by choosing nodes with a probability proportion- 
nal to the inverse of their degree in the rewiring process. 
This would account for the fact that people may be aware 
that being in contact with an highly-connected individ- 
ual is more dangerous than with someone having only 
few acquaintances. This feature could potentially hinder 
to some extent the formation of high-degree susceptible 
nodes in the network. 
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Appendix: The formalism of Gross et al. 

Until now, existing models of epidemic spreading on 
adaptive networks have been studied analytically with 
the help of low order moment-closure approximations ^51 
[^5H?7] . Here we present the formalism of Gross et al. for 
SIS dynamics on adaptive networks [53]. 

Let [Xi], [X1X2] and [X1X2X3], where X, e {S,I}, 
represent the zeroth, first and second order moments of 
the system. [Xi] correspond to the fraction of Xi nodes 
in the network, [X1X2] is the density of ^1X2 links per 
node and [X1X2X3] is the density of XiX2Ar3 triplets 
per node in the networks. The variables used by Gross et 
al. then relate to those in our formalism by the relations 
[S] = S, [I] - /, [SS] ^ Ss/2, [II] -' I1/2, [SI] = 
Sj = Is, [SSI] — Ssi, and [I SI] = Sjj. Conservation 
relations (l7| and (Is]) can now be written as [S] + [/] = 1 
and [SS] + [SI] + [II] = {k)/2. In addition to the latter 



constraints, the dynamics of the zeroth and first order 
moments of the system are described by the following 
balance equations: 



m 

dt 

dm 

dt 

d[SS] 

dt 



= P[SI]-a[I], (A.l) 

= p{[SI] + [ISI])-2a[II], (A.2) 

= {a + -i)[SI]~P[SSI]. (A.3) 



This dynamical system captures the effect of rewiring via 
the first term in (A.3 1. However, it does not yet represent 



a closed model because of the appearance of the second 
order moments in the last two equations. For this rea- 
son, the pair approximation technique is used. The latter 
consists of approximating the second order moments by 

[X1X2X3] w [X1X2] [X2X^]/[X2], which gives in this case 
[I SI] w [Slf/[S] and [SSI] w 2[SS][SI]/[S]. Together 
with this approximation, Eqs. (A.1|-(A.3) and the two 



conservation relations constitute a closed model which 
can be studied in the framework of nonlinear dynamics. 
The main difference between our improved compart- 
mental formalism and previous approaches resides in the 
level of coarse-graining in the system. The variables Ski 
and /fc/, on which our formalism is based, simply cor- 
respond to the underlying distributions of the moments 
[Xi] and [X1X2]. 
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